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In this paper, a rigorous computational method to enclose eigendecompositions of 
complex interval matrices is proposed. Each eigenpair x = (A, v) is found by solving 
a nonlinear equation of the form f{x) = via a contraction argument. The set-up of 
the method relies on the notion of radii polynomials, which provide am efficient mean 
(N of determining a domain on which the contraction mapping theorem is applicable. 
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1 Introduction 



Computing eigenvalues and eigenvectors of matrices is a central problem in many fields of 
applied sciences involving mathematical modelling. When applied to real-life phenomena, 
any model needs to consider the occurrence of diverse errors in the data, due for instance to 
lO inaccuracy of measurements or noise effects. Such uncertainty in the data can be represented 

by intervals. In the context of studying a matrix with uncertain entries, one can consider an 
^-H interval matrix A whose entries consist of intervals containing all possible errors. Although 

^-H the entries of A are intervals, many classical questions of linear algebra can be raised. For 

instance, can we demonstrate, given an interval matrix A, that for any matrix A G A, A 
^ is invertible, diagonalizable or can we enclose rigorously its eigendecomposition, that is its 

*^ set of eigenvalues and eigenvectors? In order to address explicitly one of these question, 

this paper aims at computing rigorous enclosure of eigendecompositions of interval n x n 
complex valued matrices A. More precisely, we propose a method to construct (if possible) 
a list of balls {Bi : i = 1, . . . , n} such that, for every A G A and for each i = 1, . . . , n, 
there exists Xi = (Ai,Vi) G Bi such that Av, = AjV,, where the solution Xi is unique up to 
a scaling factor of Vj. 

Before proceeding further, we hasten to mention that several different methods have been 
developed to find error bounds for computed eigenvalues and eigenvectors of standard (non 
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interval) matrices (e.g. see PQ[51|3]). Also, while the problem of computing rigorous bounds 
for the eigenvalue set of interval matrices is well studied, see for instance [H [S] and the 
references therein, a not so large literature has been produced regarding the simultaneous 
enclosure of the eigenvalues and eigenvectors of interval matrices. In this direction we refer 
to [U |6] , where different techniques have been developed to enclose simple eigenvalues and 
corresponding eigenvectors, while for double or nearly double eigenvalues a method has been 
introduced in [7 . For the rigorous enclosure of multiple or nearly multiple eigenvalues of 
general complex matrices, a significant contribution has been made by S. Rump in ^3, 8J. 

In this paper, we propose the new idea of enclosing rigorously the eigendecomposition of 
complex interval matrices by using the notion of radii polynomials, which provide a compu- 
tationally efficient way of determining a domain on which the contraction mapping theorem 
is applicable. The radii polynomials approach, which is very similar to the Krawczyk opera- 
tor approach, aims at demonstrating existence and local uniqueness of solutions of nonlinear 
equations. The approach involving the Krawczyk operator consists of applying directly the 
operator to interval vectors (in the form of small neighbourhoods around a numerical approx- 
imation) and then attempt to verify a posteriori the hypotheses of a contraction mapping 
argument [SI 110) . On the other hand, the radii polynomials are a priori conditions that 
are derived using analytic estimates, and once they are theoretically constructed, they are 
used to solve for the sets (also in the form of small neighbourhoods of a numerical solution) 
on which a Newton-like operator is a contraction. The advantage of this approach is that 
most of the estimates are done analytically and generally, and that costly interval arithmetic 
computations are postponed to the very end of the proofs. It is worth mentioning that the 
radii polynomials were originally introduced in [llj to compute equilibria of PDEs with the 
goal of minimizing the extra computational cost required to prove existence of solutions of 
infinite dimensional PDEs 12^. 

In this paper, it is demonstrated that in the context of computing rigorous enclosure of 
complex matrices, the method based on the radii polynomials is faster than the algorithm 
introduced in with the extra advantage of having local uniqueness (see Section 3.3). Also, 
it is demonstrated that an approach based on the Krawczyk operator can be significantly 
slower than the method involving the radii polynomials (see also Section 3.3). As in the 
case of the Krawczyk operator, the radii polynomials verifies existence and uniqueness of a 
zero of a nonlinear function within the inclusion interval. This also implies multiplicity 1 of 
the solution as well as non-singularity of the Jacobian at the solution. Therefore, our new 
proposed method will necessarily fail to enclose multiple or nearly multiple eigenvalues, a 
constraint that the method proposed in [5| does not have. 

The paper is organized as follows. In Section [2j we introduce the computational method, 
where we first present the method for non interval matrices. In Section [2. 1[ we demonstrate 
how to generalize the idea to rigorously enclose eigendecompositions of interval matrices. 
Finally, in Section [3j we present applications of our method. In Section |3.1[ we use the 
method to rigorously compute the Floquet exponents of a periodic orbit of the Lorenz 
system of ordinary differential equations. In Section 3.2 we study the applicability of our 
approach to matrices with interval entries of large radius. Finally, in Section 3.3 we evaluate 
the cost of our method and compare it to the cost of the algorithm of S. Rump introduced 
in [5] and to a method based on the Krawczyk operator. 



2 The computational method 

To begin with, let us fix some notation: throughout this paper we denote by IC"^" the set of 
complex matrices with interval entries, A E C"^" an n x ti complex matrix and A E IC"^" 
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an n X n interval complex matrix, meaning that any entry of A is a complex interval of the 
form 

The matrix A € C">^" ig called the center of A while ra^^j, ^'ad^^^ are called the radii of the 
real and imaginary part of Ak.j, respectively. A matrix A is said to belong to A, denoted 
A G A, ii A}^ j e Ak j for any 1 < k,j < n. Bold face letters will always denote interval 
quantities. Moreover, unless differently specified, 

• |-| is the complex absolute value and, in case of matrices M € C"^™, it acts component- 
wise, i.e. \M\i,j = |Mjj|; 

• given two real matrices Af, N any relations <, >, < etc, is assumed component-wise; 

• /„ denotes the n dimensional identity matrix, 1„ is the column vector of length n with 
all the entries equal to 1; 

• given any matrix M E ([^nxm^ ^^-le object (A/)j, stands for the n x (m — 1) matrix 
obtained by deleting the fc-th column of M. 

As already mentioned in the Introduction, given A e IC"^", the goal of this paper is to 
develop a computational method to construct a list of sets {Bi : i — 1, . . . ,n} such that, for 
any A € A and for any i = 1, . . . , n, there exists (A, v) S Bi solving Av — Av, with v unique 
up to a scaling factor. 

To simplify the exposition, we first present the method in the context of non interval 
matrices A S C"^", that is one introduces a method to enclose the solutions (A,v) of the 
equation 

Ay = Av. (1) 

As one shall see in Section |2.1[ only minor modifications are necessary for the extension to 
the interval case. 

Suppose that an approximate eigenpair of A has been computed, that is (A,v) such 
that Av « Av and let f{x) be the function / : C" — > C" that maps a point x — 

(A,Vi,V2, . . . ,Vfc_i, Vfc+i, . . . , v„) to 





Vl 




Vl 


m = A 


Vfe 


-A 


Vfc 




. v„ 




. v„ 



(2) 



where Vfc is the largest component of v. By definition, a solution x of f{x) = corresponds 
to an eigenpair (A, v) of A with the eigenvalue A given by the first component of x and the 
eigenvector v = (vi, . . . , Vfe_i, Vfc, Vfc+i, . . . , v„). Thus, the target of the incoming analysis is 
to prove the existence and to provide rigorous bounds for the zeros of f{x). Note that the 
unknowns in the equation f{x) — are A and n — 1 component of v, while the remaining 
component Vfc is a fixed parameter of the problem. Since the eigenvectors are invariant 
under rescaling, the solutions of ([T]) come in continuous families. However, fixing one of 
the component of v (in our case letting Vfe = Vfe) removes such arbitrariness and therefore 
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isolates the zeros of /. Denoting x = (A, vi,V2, . . . ,Vfe_i,Vi;_|_i, . . . , v„) and Df{x) the 
Jacobian matrix of / at one has that 



Df{x) = 



Vl 



(3) 



The problem of finding the zeros of the function f{x) is addressed by introducing an operator 
T on a Banach space whose fixed points correspond to solutions of f{x) = 0. 
Let O be the Banach space 

O = {a; e C" : ||a;||a < oo}, ||a;||n = max{|a;j|}, 



and define the operator T : — >■ by 



T{x) =x- 



(4) 



where i? is a numerical inverse of Df{x), i.e. R ■ Df{x) « /„. Since fixed points of T 
correspond to zeros of f{x), the idea is to construct a small set S C f2 such that T : B ^ B 
is a contraction, and then to apply the contraction mapping theorem to conclude about the 
existence of a unique fixed point of T in B. Note that x is an approximate zero of / and 
the operator T has been defined as a Newton-like operator around the point x, thus it is 
advantageous to test the contractibility of T on neighbourhoods of x in O. 
More precisely, denote by 

B{r) = {xen,\\x\\n<r} 

the close ball of radius r around the origin and let Bx{r) = x + B{r) the ball with the 
same radius and centered at x. Treating r as a variable, we choose the balls Bx{r) as the 
candidate sets where to check if T is a contraction. The question whether T is a contracting 
map will be formulated in terms of the verification of a set of computable conditions, called 
radii polynomials, whose construction is based on verifying efficiently the hypothesis of the 
next result. 

Theorem 2.1. Suppose that Y, Z{r) G fl are such that 

\T{x) -x\<Y, sup \DT{x + b)c\ < Z{r), (5) 

b,ceB{r) 

and satisfy 

\\Y + Z{r)\\n<r. 
Then there exists a unique x € B^ir) such that T{x) = x. 

Proof. The mean value theorem applied component- wise to T implies that for any x,y G 
Bx{r) and for any A; = 1, . . . , n, 

Tk{x) - Tk{y) = DTk{z){x - y), for z e {tx + {1 - t)y : t G [0, 1]} C Bs{r). 



Then, 



\Tk{x)-Tk{y)\ = 



DTk{z) 



r{x - y) 



\\x-y\\n 



-\\^-y\\n<^\\x-y\\n. 
r r 



(6) 
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Choosing y — x and using the triangular inequahty, one has that 

\Tk{x) -Xk\< \Tk{x) - Tk[x)\ + \Tk{x) -Xk\<Yk + Zk{r) < r, 

where the last inequality follows from the fact that Y, Z g M", and thus ||F + 2'(r)||n = 
maxi{Y, + Zi{r)}. That proves that T{Bs{r)) C Bg(r). From (|^, it follows that 

\\T{x) - T{y)\\n = max{|r,(x) ~ n{v)\} < - y\\n. 

k r 

Since ||Z(r)||s7 < r, T is a contraction on Bx{r). Thus, from the contraction mapping 
theorem, the exists a unique fixed point of T in Bx{r). 

□ 

Definition 2.2. Given the vectors Y, Z{r) S satisfying ([5]), we define the radii polynomials 
Pk{r), A; = 1, . . . ,n by 

Pkir) = iY + Z{r))k-r. 

The following result holds. 

Lemma 2.3. Consider the radii polynomials Pkif), k = l,...,n introduced in Defini- 
tion\2.S\ Then for any r > such that 



Pk{r) < 0, for all fc = 1, . . . , n, 
there exists a unique x e Bx{r) such that f{x) = 0. 

Proof. Suppose f > is such that Pkir) < 0, for all k — 1, . . . ,n. For ^ all the entries 



of Y, Z{r) are real positive, thus maxi{{Y + Z{r))i} = \ \Y + Z{r)\\ii < r. From Theorem 2.1 
there exists a unique x G Bx{r) such that T{x) — x and therefore f{x) — 0. 

□ 

We proceed with the explicit construction of the bounds Y and Z. Since T{x) — x = —Rf(x), 
let 

Y \Rfix)\. (7) 
The bound Z{r) satisfying ([s]) is constructed as a polynomial in r. First rewrite DT{x + b)c 

DT{x + b)c = {In - R ■ Df{x))c + R[{Df{x) - Df{x + b))c] 



as 



so that 
Define 
where 



\DT{x + b)c\ <\iIn-R- Dfix))c\ + \R[{Df{x) - Df{x + b))c]\. (8) 
Z{r) = rZ^+r^Z^, (9) 



Zo = \In~R- Df{x)\tn. Zi = 2|i?|l„, (10) 
where 1„ is the same as 1„ with a zero, instead of one, in the fc-th component. 
Lemma 2.4. Consider the polynomial vector Z(r) defined by (|9|. Then 

Z{r)> sup \DT{x + b)c\. 

b,ceB{r) 

Proof. From (Isl), the statement follows by proving that 
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i) sup \{Ir,-R-Df{x))c\<rZo 

ci£B{r) 

ii) sup \R[{Df{x)~Df{x + b))c]\<r^Zi. 

b,ceB(r) 

Since c S B{r), \c\ < rln, it follows that 

sup \{In-R-Df{x))c\<rZo . 

That proves i). For any b = (&a, • ■ ■ , ^fc-i, &fe, • • • , 





&1 















(bxinh 




bk+i 






bn 





Note that the fc-th row of the above matrix is null. Since \bi\ < r, we have \ {Df{x) — Df{x- 
b))c\ < 2rH„ and therefore supft_^gB(r) \R[iDf{x) - Dfix + b))c]\ < 2r^\R\i„ = r^Zi. 



□ 



In summary, given an approximate eigenpair (A,v), the method consists of computing 
rigorously the bounds Y and Z{r) given by ([7| and ([9]), and then to check whether there 
exists an interval I where all the polynomials pk (r) are negative. If I we select r — inf X 
and we conclude that f{x) = has a unique solution within the ball B^{r). In practice, we 



get the existence of an eigenpair (A,v) of A, with |A — A| < r, 



< r, for j ^ k and 



Vfe = Vfc. To prove the existence of a second eigenpair of A, it is necessary to provide an 
approximate solution (A, v), different from the previous one, and to repeat the computation. 



2.1 Extension to the interval case 

Besides few modifications necessary to deal with interval quantities, the procedure to com- 
pute rigorously bounds for the eigendecomposition of an interval matrix A g IC"^" is 
basically the same as for the scalar case. However, a fundamental difference consists in the 
fact that all the computations are done in the interval arithmetic regime [1], in which any 
of the basic operations o G {+,—,-,/} is extended to the interval case in order to satisfy 
the general assumption 

VPeP VQeQ, PoQePoQ. (il) 

Given A an interval complex valued matrix, we now address the problem we stated 
at the beginning of the paper, that is how to construct sets {-Bi}r=i so that each A £ A 
admits one and only one eigenpair (A,v) in any of the Bi's. Recall that A is the center 
of the interval matrix A. We first compute (A, v) an approximate eigenpair of A and, as 
before, define x = (A,vi,V2, . . . , Vfe_i, v^+i, . . . ,v„) where the missing component is chosen 
so that v/j = maxj{vj}. Then, replacing the scalar matrix A in (12) by the interval matrix 
A, the function f{x) and the Jacobian matrix Df{x) defined in (2| and ([3| are replaced 
respectively by / : C" — > IC" and by an interval matrix Df{x) that represents a linear 
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operator from C" to IC". We choose i? to be a numerical inverse of Df{x), the center of 
Df{x), and we proceed to the definition of the operator T{x) — x — Rf{x) and to the 
bounds Y and Z(r), as done before with the boldface quantities in place of the previous 
one. Clearly the quantities on the left hand side of relations ([s]) are now intervals, thus we 
define component-wise Y, Zi, Z2 as the suprcmum of the intervals appearing on the right 
hand sides of ([t]) and ( [lO| . That will yield the uniform bounds 

\T{x)-x\<Y, sup \DT{x + b)c\<Z{r), 

b,c£B{r) 

where DT{x) =I„-R- Df{x). 

Suppose that r > is such that the radii polynomials Pk{r) < for all k. Then interval 
arithmetic insures that for all A E A, there exists a unique (A,v) such that |A — A| < r, 
|vj — Vjl < r, Vfe — Vfc, and Av = Av. In other words, r is a uniform bound in A for the 

existence of an eigenpair of any A E A. Indeed, having fixed (A,v) and R ~ (Df{x)) ^, 
for any A £ A define /a (2;) and DfA{x) as in pi and (|3|, and the fixed point operator 



Ta{x) = X — RfA{x). The fundamental inclusion (11 1 implies that fA{x) € f{x), DfA{x) £ 
Df{x) and Ta{x) € T{x), for any A G A and x G C". Thus, as A varies in A, the bounds 
([5]), with Ta in place of T, are satisfied for the same Y, Z, r proving the existence of a fixed 
point in Bx{r) for any Ta and consequently an eigenpair for any A Cz A. 



3 Results 

In this section we report some computational results. All the computations have been done 
in Matlab supported by the package Intlab [13, where the interval arithmetic routines have 
been implemented. The approximate eigenpairs (A, v) of A have been computed running 
the standard eig.m function in Matlab. In order to avoid rounding error and to obtain 
rigorous results, we emphasize that the computational algorithm treats any matrix as an 
interval matrix. Thus, even if one wishes to deal with a scalar matrix A, the method 
first constructs a (narrow) interval matrix around A and perform all the computation with 
interval arithmetics. 

3.1 Example 1: rigorous computations of Floquet exponents 

The first example concerns the rigorous enclosure of the Floquet exponents and related 
eigenvectors associated to a periodic orbit of the Lorenz system 

til = Cf{u2 - Ml) 

U2 = pUl -U2~ U1U3 . (12) 
"3 = U1U2 - I3U3 

For a choice of the parameters a — 10, /3 — 8/3, p = 20.8815 the existence of a periodic orbit 
7(t) in a neighborhood of an approximate solution has been proved [H]. It is known that 
the stability character of the orbit ^{t) is encoded by the Floquet exponents, which are the 
eigenvalues of a particular real matrix, denoted by A, resulting by integrating the linearized 
system around 7(i). Without going into the details (we refer to [TS] for an exhaustive 
explanation), we only mention that one of the Floquet exponents is zero, due to the time 
shift, while the number of the other eigenvalues of A with negative (positive) real part 
gives the dimension of the stable (unstable) manifold. Moreover the associated eigenvectors 
provide the directions tangent to the invariant manifolds on 7(t) (tangent bundles). 
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By means of a computational method based on the radii polynomials, in |15| we proved 
that the matrix A associated to the solution j(t) lies within the interval matrix A centered 
at 

-10.55360193 5.33379647 -5.24740415 
A= 0.31403414 2.33062549 -3.32865541 
-7.49045333 5.01386821 -5.44369022 

with radius rad = 9.66146973 • 10^*", meaning that each entry A{i,j) consists of the interval 
[A{i,j) — rad,A{i,j) + rad]. Following the computational method discussed in Section [i] 
we compute the enclosure of the eigenpairs of A: it results that any A G A admits three 
eigenpairs (Aj, v^), i = 1, 2, 3 each one lying in the ball of radius r.i around the approximate 
values (Ai,Vi) given in Table [ij 





i = 1 


i = 2 


i = 3 




2.774764083439355 • 10^^ 


3.567796353801448 • 10"^ 


3.649406638638561 • 10"^ 


Ai 


-13.962049368058944 

r -0.777880129851759 ' 

-0.111361799590875 
[ -0.618466695282529 


-9.363255645359691 • lO^^^ 

r 0.121332037790074 ' 

0.807699961688802 
[ 0.576974270218020 


0.295382701392333 

r 0.156626754180928 " 

0.835985628946308 
[ 0.525924032603563 



Table 1: Data associated to the rigorous computation of the Floquet exponents of the 
periodic orbit 7. 



We remark that in the general situation the genuine solution (A, v) of the eigenproblem 
is proved to exist in a complex neighborhood of the approximate solution (A, v). Therefore, 
even if one or both A and v are real vectors, the same can not be concluded for A or v. 
However, if the matrix A and the approximate solution A and v are real and the computation 
is successful, then the genuine solution so obtained by solving the radii polynomials is also 
real. Indeed, suppose the contrary, that is the exact solution A and v are complex. Since 
A is real, the complex conjugate couple (C(A),C(v))) is also a solution of the eigenproblem, 
AC{v) = C(A)C(v). But both the solutions (A,v) and (C(A),C(v))) belong to the same ball 



in 17 around x and this violates the uniqueness result stated in Lemma 2.3 The same 
argument extends in the case of interval matrices. 

Coming back to the previous example, we conclude that the Floquet multipliers of "f{t) 
are real, one is negative and one is positive. Note that the zero Floquet multiplier is indeed 
contained in the ball around (A2,V2). 



3.2 Example 2: matrices with interval entries of large radius 

In the next sample we compute the eigendecomposition of an interval matrix A constructed 
as follows: consider the complex number Aq = and Xj — e*~-', j — 1, ... ,5 and define D as 
the diagonal matrix with entries A^, A: = 0, . . . , 5. Let A — XDX^^, for a random matrix X 
with values in the complex square [— 1, 1] + i[— 1, 1] and finally let A be the interval complex 
matrix centered at A with component-wise radius rad both in the real and imaginary part. 
For different values of rad we compute the enclosure of the eigenvalues of A: given Afc the 
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Figure 1: Balls in the complex plane enclosing the 6 eigenvalues of any A £ A for rad = 1.3 ■ 10 



approximate eigenvalues of A, 



Ai - 

A2 = 
A3 = 
A4 = 
A5 = 



0.00000 
0.30901 
-0.80901 
-0.80901 
0.30901 
1.00000 



O.OOOOOi 
0.95105i 
0.58778i 
0.58778i 
0.95105i 
O.OOOOOi 



denote by r^, fc = 0, . . . 5 the radius of the ball in the complex plane centered at inside 
which, for any A G A, a unique eigenvalues of A has been proved to exist. The results 
are presented in Table [2j where the different values of rad are given in the first column, 
while each row collects the result for the radii r^. For values of rad < 1.3 • 10~^ the method 
succeeded in computing the enclosure of the entire eigendecomposition of A, while for larger 
values of rad the method starts failing in enclosing some of the eigenpairs, up to the value 
rad = 3.2 • 10~^, where no computation is successful. Figure [l] shows the radii of the disks 
in the complex plane enclosing the six eigenvalues for the case rad = 1.3 • 10^'^. 



3.3 Example 3: evaluating the performance 

In this last section we compare the performance of our new method, that we denote here by 
radiipol, with two different algorithms developed by S. Rump. The first one, here denoted 
by verifyeig, has been introduced in with the primary goal of computing enclosures of 
multiple of nearly multiple eigenvalues (and related eigenvectors) of interval matrices. It con- 
sists of a verification method and provides rigorous bounds around an approximate solution 
within which the existence, but not uniqueness, of the exact solution of the eigenproblem is 
proved. The second one, here denoted by verifynlss, is based on a Krawczyk operator [9l[T0] 
and is a general routine to rigorously compute well separated zeros of nonlinear functions. 
In fact, in the code verifyeig.m (available in the library Intlab |13j). where the method 
verifyeig has been implemented, the author suggested to use verifynlss to compute simple 
and well separated eigenpairs. This method is implemented in the code verif ynlss .m in 
the library Intlab [T5] . 
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rad 


fo 


ri 


r2 


rs 


r4 


r5 


1 • 10"^ 


2.0 • 10-4 


1.5 • 10-4 


1.6- 10-4 


1.3- 10-4 


1.8 • 10-4 


1.6 • 10-4 


1 • 10-4 


0.0021 


0.0016 


0.0018 


0.0014 


0.0019 


0.0017 


1 • 10-3 


0.0281 


0.0181 


0.0202 


0.0149 


0.0219 


0.0197 


1.3 • 10-3 


0.0492 


0.0247 


0.0284 


0.0201 


0.0306 


0.0276 


1.9 • 10-3 




0.0416 


0.0563 


0.0323 


0.0515 


0.0581 


2 • 10-3 




0.0452 




0.0346 


0.0678 


0.0591 


2.5 • 10-3 




0.0759 




0.0485 






3.1 • 10-3 








0.0828 






3.2 • 10-3 















Table 2: Enclosures of the eigenvalues of a complex interval matrix A, as the radius rad of 
the entries of A increases. 



We summarize the obtained results in Table |3] and Table [4] In Table [3j the computa- 
tional time necessary to compute the entire eigendecomposition of a test matrix (measured 
in seconds on a 2.4 GHz computer using the tic-toe Matlab function) is reported, inde- 
pendently whether the methods were successful or not. Table [4] presents the output of the 
algorithms: for each computation of the eigendecomposition, it provides the average of the 
radius of the balls enclosing the exact eigenpairs. The entry — means that the method fails 
in the enclosure of at least one of the eigenpair. 

For both experiments the test matrices A have been constructed as in the previous 
section: given TV we define D e £N+i,n+i ^ diagonal matrix with entries given by A'' 
equispaced values on the unit circle in the complex plane and 0, i.e. diag{D) = [0, e^^^],j = 
1, . . . ,N. Then let A = XDX~^, where X is a complex random matrix with entries in the 
complex square [— 1, 1] 1, 1] and finally define A as the interval complex matrix centered 
in A and of radius rad. 

From the results of Table |3j we conclude that our new proposed method radiipol is much 
faster than the method verifynlss based on the Krawczyk operator, while also proving 
existence and uniqueness of the solutions. Moreover, our algorithm is faster than verifyeig, 
especially for small N, and the computational time of the two methods is almost equivalent 
for large N. This achievement have been possible thanks to the analytical estimates (in the 
form of the radii polynomials) introduced in the method, that allow minimizing the number 
of computations done with interval quantities. Hence, the radiipol method provides a 
computationally efficient way to compute eigendecompositions of complex interval matrices. 

The results presented in Table |4] confirm that the new approach radiipol is satisfactory 
also from the point of view of the accuracy of the results. Indeed, while the algorithm 
verifynlss fails quite soon as TV and rad increase (it fails for rad = and for all N > 15), 
the new algorithm is successful also for large entries of A. Moreover, as one can see in 
Table |4] the performance of radiipol is very close to the performance of the algorithm 
verifyeig which, we underline this one last time, does not prove uniqueness of the solution. 
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Computational time (s) 



N 


radiipol 


verifyeig 


verifynlss 


5 


0.0148 


0.0372 


0.3423 


10 


0.0282 


0.0650 


0.6432 


50 


0.2481 


0.4170 


7.4476* 


100 


1.4580 


1.8814 


66.119* 


200 


14.512 


16.739 


792.76* 


500 


377.3 


423.6 





Table 3: Comparison of the running time necessary to enclose the entire eigendecomposition of 
A using radiipol and the two other methods verifyeig and verifynlss. A is a complex interval 
matrix of dimension + 1 and rad = 10"^^. The entries of the form (•)* correspond to unsuccessful 
computations. The figure on the right depicts the ratio between the computational time of the 
algorithms verifyeig and radiipol. 



Average of the radius of the disks enclosing the eigenvalues 





N=5 


N=10 


rad 




10-^" 


10- 


4 


10-"' 


10-5 


lO-'' 10-'' 


radiipol 




2.7640"^^ 


0.0019 


3.25T0- 


4.61T0- 


4 




verifyeig 


4.69 TO-15 


2.07T0~9 


0.0016 


2.08T0- 


'9 3.02-10- 


" 0.0049 


verifynlss 


6.26 -10-^ 


N=50 








N=100 




N=150 


rad 


io-^» 


io-« 






10-"' 


io-« 


10-' 


10-iu 


radiipol 


2.69T0-'^ 


4.94-10-5 




9.02-10-'^ 






1.31-10-*^ 


verifyeig 


5.59T0-8 


9.45-10-s 




1.3140-'^ 


2.07-10-5 




1.64-10-^ 


verifynlss 



















Table 4: Comparison of the accuracy of the three methods as the dimension A'^ and the radius rad 
of the test matrix A change. The entry — means that the method fails in the enclosure of at least 
one of the eigenpair. 
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